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Computer Algebra (CA) 



In order to get an idea what facilities CA provides it might be sensible to ask what we 
expect from CA. I suggest that 

CA is supposed to solve or simplify any algebraic problem that we can 
state explicitly. 

So, CA has to provide methods to simplify algebraic expression and to solve equation 
systems. This, of course, is not enough because it also has enable us to input our problem 
and has to generate some output of the result (also graphical, plots, LaTeX, etc.). I 
want to differentiate the needs for the input a little bit more: On the one hand, CA has 
to offer some language or formalism that includes many operators to form expressions 
(+, *, sin, f, ^2, etc.) and on the other it should provide methods to make definitions, 
declarations, substitutions, and assumptions. Thus, we split the facilities we expect from 
CA into four groups: 

A operators to form expressions 

B methods to make definitions, declarations, substitution, and assumptions 

C methods to simplify and solve 

D methods to generate output 

The reason why I suggest this splitting of the commands is the following: It seems that 
the discussion of group A dominates most handbooks on CA only because it is by far the 
largest group of commands and operators. I find though, that the commands in group B 
and C are more elementary and important to understand whereas for the commands in 
group A one can well refer to the online user's manual. The commands in group D are 
also important for efficient work with CA. This is why these lectures try to emphasize the 
commands in groups B, C, and D a little bit more than usual. 
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Chapter 1 



Reduce I 



1.1 The simplification principle 



As stated, I think the simplification and solving methods of a CA system to be most 
elementary and important. Hence we start to discuss these in the case of Reduce. 

Reduce is an input-output machine. One could say that all that Reduce does is to 
reformulate your input expression obeying certain rules! One of these rules is, 
e.g., to execute all operators and commands in the expression. For example, Reduce 
reformulates 1 + 1 by executing the + operator and answering 2. Fortunately, the rules 
for reformulation are such that - usually - this reformulation means a simplification. 
These rules can be influenced by the user: Either by switching on and off the rule switches 
or by introducing new rules. Most important, Reduce offers the rule switches collected in 
table 1.1. 

As one can see, by switching on the exp-rule, Reduce reformulates and thus simplifies the 
expression (x + l)(x — 1) as x 2 — 1. This principle is quite different to other CA systems. 
The combination of such rules can be very powerful. The rule switches turned on at 
default are: allfac, exp, mcd, 1cm, ratpri, pri, nat. This means that, at default, 
for any expression (subexpression, component, etc. ), Reduce expands the multiplication 
of two larger terms, interprets all divisions as rationals, tries to divide common factors in 
numerators and denominators, and finally tries to factorize simple factors. This is already 
quite a powerful simplification scheme. 
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switch 


description 


e.g. 


allf ac 


factorize simple factors 


2x + 2 -> 2(x + l) 


div 


divide by the denominator 


(x 2 + 2)/x -> x + 2/x 


exp 


expand all expressions 


(x + l)(x-l) ^x 2 -1 


mcd 


make (common) denominators 


x^ 1 — > 1/x 


1cm 


cancel least common multiples 




gcd 


cancel greatest common divisor 




rat 


display as polynomial in factor 


z+l^l + x- 1 


ratpri 


display rationals as fraction 


x 1 

' X 


pri 


dominates allfac, div, rat, revpri 




revpri 


display polynomials in opposite order 


X 2 + X+ l^l+X + X 2 


rounded 


calulate with floats 


1/3 -> 0.333333333333 


complex 


simplify complex expressions 


l/i -> -i 


nero 


don't display zero results 


-> 


nat 


display in Reduce input format 


sd x **2/3 


msg 


suppress messages 




fort 


display in Fortran format 




tex 


display in TeX format 





Table 1.1: Switches for Reduce's reformulation rules. Those marked with * are turned on 
at default. 



1.2 The interface 



Starting Reduce on Unix it displays about 

Loading image file: 

/vol/sc/lib/reduce3.6_patch980830/reduce.img REDUCE 3.6, 15-Jul-95, 
patched to 30 Aug 98 . . . 

1: 

The first three lines are an artifact of Reduce being implemented in LISP: The program 
code (image) of Reduce is passed to the LISP interpretor. The prompt 1 : expects a 
command as input. Reduce ignores all upper case letters! A command consists out 
of a statement and a terminator. The terminator decides whether the Reduce's response is 
displayed (; terminator) or not ($ terminator). The following lines are command lines only, 
the response of Reduce is (sometimes) indicated after the comments sign °/„. All examples 
in this chapter are collected in the file redl at [5]. 

7, file "redl" 

off factor, exp, mcd, ratpri$ 

(x+l)**3/3 - x; '/,-> 1/3* (x + 1)**3 - x 

on mcd$ ws; '/,-> ((x + 1)**3 - 3*x)/3 

on ratpri$ ws; %-> 
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on exp$ ws ; 

factor x$on rat$ ws ; 
remfac x$off exp$ 



'/.-> (x**3 + 3*x**2 + l)/3 
%-> x**3/3 + x**2 + 1/3 



z+a; 

order z,a$ ws; 
order nil$ ws; 



%-> a + z 
°/.-> z + a 
'/.-> a + z 



1.3 User definitions 



1.3.1 Names and assignments 

off mcd$ 

my_namel23 := (x**2-l)/(x+l) ; '/.-> my_namel23 := (x + 1) ** (-1) * (x**2 - 1) 
on mcd$ my_namel23 ; %-> x - 1 

clear my_namel23; 

A string can be used as a name (or identifier) for variables, procedures, etc. The first 
character of a name has to be alphabetic, whereas others can be numbers or e.g. 'all' 
and 'ric_scalarO' are proper names. A name should not coincide with a reserved variable 
name (e i infinity nil pi t) or the command names in table 1.2. 

As long as a name is not assigned to some (symbolic) value it is considered to be a clear 
name. The assignment operator 

name: =expr; 

assigns the (reformulated) expression on the rhs to the name on the lhs. To withdraw an 
assignment one uses the command 

clear name; 

This command also deallocates the memory resource that is associated with the name. 

For understanding Reduce, it is important to distinguish between the expres- 
sion assigned to a name and the value that Reduce displays when evaluating the 
name. The last of which was reformulated by Reduce with current rules and definitions: 

a:=b$ a; %-> b 

b:=l$ a; %-> 1 

a:=a$b:=2$ a; %-> 1 

clear a,b; 

In the first line, a is assigned to b and is, of course, also evaluated to be b. In the second 
line, a is still assigned to b but is evaluated to be 1 (because b is assigned to 1). In the 
third line, a is assigned to 1 (because the rhs a is evaluated to be 1) and an thus it is also 
evaluated to be 1 even after b is assigned to 2. We see that the reassignment (a:=a$) is not 
trivial at all! Analogously, the following example shows the importance of the reassignment 
if rules have changed: 
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off mcd$ a: = (x**2-l) / (x-1) ; '/.-> a := (x**2 - l)*(x - 1)**(-1) 

on mcd$ a := a$ a; '/,-> x + 1 

off mcd$ a; '/,-> x + 1 

Here, in the end, a is assigned to x+1. Without the reassignment a:=a$, a would still be 
assigned to (x**2 - l)*(x - l)**(-l) and also evaluated to be such in the last line. 

The reassignment is one of the most important tools to apply new definitions 
and to use the simplification methods of Reduce. 



1.3.2 Aliases 

To define aliases for any expression, one can use the define command. The syntax is 

define alias= expression; 

Note that the rhs expression is evaluated before it is assigned to the alias. However: An 
alias can never ever be changed or unassigned during a Reduce session! If an 

alias appears in an expression it is replaced by its associated value before Reduce applies 
any other rules. Example: 

define isprime = primep; 

for x:=0:20 do if isprime (x) then write x," is a prime!"; 

'/. throws 2 3 5 7 11 13 17 19 

define is_assigned_to = :=; 

x is_assigned_to 5; '/,-> x := 5 



1.3.3 Substitutions 



Actually, we already noticed how to substitute identities into expressions: Define the iden- 
tity with an assignment statement or introduce a rule for this identity. Then reevaluate 
the expression. Say, e.g., we want to substitute x = into f(x) = cos(x), then we could 
write one of the following three lines: 



cos(x)$ x:=0$ f ; clear x$ '/,-> 1 

cos(x)$ let { x=>0 }$ f ; clear x$ '/,-> 1 

cos(x)$ f where { x=>0 }; '/,-> 1 

cos(x)$ sub(x=0,f); '/,-> 1 



The last possibility is the most beautiful, because one does not have to clear the assignment 
x : =0 or the rule x=>0 hereafter. The where command only introduces a local rule. After 
these substitutions, f is still assigned to cos(x) but was only locally evaluated to be 1. 
Reduce also offers the match command to substitute expressions for polynomials instead 
of names. 
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1.3.4 Rules 

In the upper example we already introduced the let command that introduces new rules. 
In general, a rule has syntax 

expr => expr [when boolean] 

and the let command accepts a list of rules as parameter. Rules can be of explicit nature 
where any literal occurrence of the lhs expression is replaced by the rhs expression. Or 
they can by of parametric nature where the lhs expression includes formal names (like 
parameters) that represent any expressions. These formal names are marked with a twiddle 
such as ~x: 

sin(x)**2 + cos(x)**2; %-> cos(x)**2 + sin(x)**2 

trig_rul.es : ={ 

cos (~x) *cos (~y) => (1/2) * (cos (x+y)+cos (x-y) ) , 

sin(~x)**2 + cos(~x)**2 => 1 }$ 
let trig_rules$ 

sin(x+l)**2 + cos(x+l)**2; %-> 1 

cos(pi/3-l) * cos(l); %-> (2*cos((pi - 6)/3) + l)/4 

showrules cos; 7, shows ALL rules associated with cos! 

clearrules trig_rules; 

As already mentioned, the where command applies rules only one the preceding expression. 
The command 

showrules name; 

displays all rules associated with the name. This gives an insight in how Reduce defines 
functionals: 

1: showrules log; 
-Clog(l) => 0, 
log(e) => 1, 
~x 

log(e ) => x, 

1 

df (log(~x) ,~x) => — } 
x 

1.3.5 Operators 

An operator is any name that accepts parameters. In Reduce, one can declare a name to 
be an operator without at all specifying how many parameters the operator expects and 
what the functionality of the operator is! This leaves the user large freedom to handle 
and use operators. E.g., one can introduce indexed coordinates x l by declaring x to be an 
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operator with the index as only parameter and without any further properties. In fact, 
most other, more specialized declarations (like matrix, procedure, etc.) can be thought of 
special kinds of operator declarations. And also functions (such as sin, log, etc.) can be 
introduced as operators with certain properties. The properties of an operator are declared 
by specifying the rules how Reduce has to handle them. Example: 



operator eps$ 
antisymmetric eps$ 
let { eps (0,1,2,3) => 1 }$ 
eps(3, 1,2,0) ; eps(2, 0,1,1) ; 
showrules eps; 



% declares a new operator 

% declares eps to be totally antisymmetric 
7, defines the properties of this operator 
7.-> -1 and 

•/.-> {eps (3, 2, 1,0) => 1}$ 



clear my_sin$ operator my_sin$ 

my_sin_rules : ={ 

my_sin(~x) => - my_sin(-x) 
my_sin(~x) => my_sin(x - 2*pi) 
my_sin(~x) => - my_sin(x - pi) 
my_sin(~x) => my_sin(pi - x) 
my_sin(0) => 0, 
my_sin(pi/6) => 1/2, 
my_sin(pi/4) => l/2*sqrt(2), 
my_sin(pi/3) => l/2*sqrt(3), 
my_sin(pi/2) => 1 

>$ 

let my_sin_rules$ 

my_sin(-124*pi/6) ; '/.-> 
clearrules my_sin_rules ; 



when numberp(x/pi) and x/pi<0, 
when numberp(x/pi) and x/pi>2, 
when numberp(x/pi) and x/pi>l, 
when numberp(x/pi) and x/pi>l/2, 



sqrt(3)/2 



1.4 Solving 

solve(x**3-3*x**2-61*x+63,x) ; '/„-> {x=9 ,x=l ,x=-7> which are the nulls of this polynomial 
solve({x+y-9,-2*y**2+3*x},{x,y}) ; '/.-> { {y= - 9/2,x=27/2> , {y=3,x=6> } 

Reduce tries to solve algebraic equation systems with 

solve ( exprJist, name-list) ; 

Without the switch multiplicities Reduce display multiple solutions only once. 

Reduce has no build-in routines to solve differnential equations. However, there exist pack- 
ages for this problem. One should especially mention the CATHODE 2 project (Computer 
Algebra Tools for Handling Ordinary Differential Equations) which develops routines to 
hande differential equation systems also for Reduce. The members of this group are M. 
MacCallum (London), V. Fairen (Madrid), E. Tournier (Grenoble), Th. Mulders (Zurich), 
M. van der Puf (Groningen), F. Schwarz (Bonn), L. Brenig (Brusssel). You find their 
packages at [10], especially Crack and ODEsolve. 
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1.5 Commands and references 



Most of the commands Reduce offers are to formulate expressions. I tried to collect most 
in table 1.2. The point of this table is to give you an idea what commands there exist 
and to briefly describe the syntax. One should read through this table once. For detailed 
explanations we refer to the user's manual which can be accessed in the internet at [6]. I 
prepared more online references at [5]. 

for 1:=1:50 sum 1; '/,-> 1275 

for 1: =1:100 product 1; 7, gives the huge number 100! 

int(sin(x)**2,x) ; '/,-> x/2 + (-cos (x) *sin(x) ) /2 

depend (f,x)$ df(f,x); int(f,x); %-> df(f,x) int(f.x) 
nodepend (f ,x) $df (f ,x) ; int (f ,x) ; %-> x*f 

matrix m(3,2)$ 7o declare matrix 
m:=mat((l,2,3) , (2,1,2) , (3,2,1))$ 7. construct matrix 

det m; 7d~>8 determinante 

1/m; 7. gives the inverse 

matrix v(3,l)$ % declare a vector 
v:=mat((l),(2),(3))$ 

1/m * v; %->( 1,0,0) matrix multiplication 

procedure facl(k)$ 

for l:=l:k product 1$ 7« defines a new factorial procedure: 

facl(16); %-> 20922789888000 which is 16! 

operator fac2$ 

let { fac2(l) => 1 , 

fac2(~n) => n*fac2(n-l) when numberp n and n > 1}$ 
fac2(16); '/.-> 20922789888000 

write "next ""prime > 500"" is " ,nextprime(500)$ '/.-> next "prime > 500" is 503 

showtime$ 7« milli seconds since last showtime 

in "my_file"$ °L executes all lines in "my_file" as input 

quit$ 7» exits Reduce 



1.6 Exercises 

1. Get used to Reduce's interactive interface. Try control commands like showtime 
write quit input ws clear °/„. 

2. Turn all switches (for the reformulation rules) off: 

off allf ac,div,exp,mcd,lcd,ratpri ,nat$ 

Now type in term (x**2 - 1) / ( (x+1) * (x-1) ) ; Which switches do you have to turn 
on additionally to get the following reformulations by Reduce: 
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control 

assignment 
unassignment 
substitutions 
rules 

elementary 
logical 

selection 
functional 

calculus 

matrix 

lists 

loops 

conditions 
groups 

ordering 

solving 

declarations 

properties 



quit showtime pause cont; ws input prompt-number; 
write exprs; in out shut loacLpackage "filename"; 
rederr "message"; °/„ comment 

name := expr; set (name, expr) ; define name:=expr; 
clear names; 

let ruleJist; sub match (eqnJist, expr) ; expr where ruleJist; 
expr => expr; let clearrules ruleJist; showrules name; 
expr +-*/** expr; 

expr = neq >>=<=< expr; ordp (names) ; 
boolean and or boolean; not boolean; 

numberp fixp evenp primep (expr); freeof ( expr, name) ; 

rhs lhs eqn; num den expr; coeff (expr, name); 
coeffn ( expr, name, degree) ; 

cos sin tan esc sec cot a" ~h a~h atan2 abs sqrt exp In log 
logb loglO hypot factorial fix ceiling floor nextprime conj 
impart repart random round sign (expr); max min (exprJist) ; 

int df (expr , names) ; depend nodepend (name, names) ; 

mat (components) ; tp trace det rank mateigen (matrix) ; 

{elements}; first second third rest reverse length (list); 
append (list, list); expr . cons list; 

for (see section 2.3.2); while boolean do statement; 
repeat statement until boolean; 

if boolean then statement [else statement] ; 

begin scalar names$ statements return expr$ end; 
<< statements »; 

order korder factor names; 

solve ( exprJist, name Jist) ; 

operator matrix names; array name( dimensions) ; 
procedure name (parameter '..names) >$ statement; 

noncom symmetric antisymmetric even odd linear operators; 



Table 1.2: Standard commands in Reduce. Explanation: Each group of commands which 
are only separated by blanks have the same syntax. These groups are separated by ' ; ' and 
the syntax in specified only once for all commands in one group (think of the blanks as 
'or'). The plural of italic words means that there can be many of them separated by 
commas. Terms in brackets '[ ]' are optional. The twidle ~ means one of sin, cos, etc. 
ea^r^expression, egn=equation (which has syntax expr=expr) 
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a) (x**2 - l)*(x + l)**(-l)*(x - 1)**(-1)$ 

b) (x**2 - l)**(-l)*x**2 - (x**2 - 1)**(-1)$ 

c) (x**2 - l)**(-l)*(x**2 - 1)$ 

d) 1$ 

We see that mcd is absolutely neccessary to calculate with rationals. The exp switch 
is important when it is not obvious that terms cancel in large expressions. Similary, 
the gcd switch improves the perfromance in canceling a common divisor of numerator 
and denominator (which mcd does also in easier cases). 

3. What is the value of 

a(a + 2) + c(c - 2) - 2ac (1.1) 
for a — c = 7? 

4. Get used to the user's manual at 

http: //www. uni-koeln. de/REDUCE/3 . 6/doc/reduce/reduce . html and to the file 
execution with the in command (See also section 2.1). Find out the syntax of the 
commands limit factorial int df for. Thereby get answers for 

lim^M, 15!, ff(x)dx, |exp(/(s)) /'(*) dz , 

i=l 

Generate the list {2°, 2 1 , .., 2 10 } (use for and collect). 

5. Does Reduce confirm 
oz\ ( oz 



dx J ' \ dy / ' 
Teach Reduce the appropriate rule for such problems! 

6. Generate a list of all prime numbers smaller than 1000 and being a divisor of 606353. 

7. Write a procedure that counts the zeros at the end of the number k\ for any positive 
integer k. 

8. Develop a procedure to determine the nth order Taylor expansion of some arbitrary 
function f(x) at one arbitrary point Test the procedure by considering e x 
and sinx. Confirm e lz = cos z + i sin z up to 5th order. 

9. Prove by complete induction that 
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Write a procedure that return the characteristic polynomial and the eigen values of 
a square matrix. (Check with but don't use mateigen.) Implement operators that 
return the scalar and vector product, a ■ b and a x b, for two vectors a, b. Use arrays 
for these problems if you want indices to run from 0. 



Chapter 2 



Reduce II 



2.1 Running files with Reduce and packages 

For most purposes it is convenient to write all commands to be executed by Reduce into a 
separate file edited with your favorite editor. The in command causes Reduces to run this 
file as if all the lines where typed in during an interactive session. Use in "filename" ; if 
you want Reduce to display each line of the file and in "filename"$ if not. The last line 
of such a file should be end$. Other useful command in this context are pause, cont, 
demo. 

On Unix systems one can also use the following command line to pipe your Reduce file 
sample. rei into reduce and collect all output in the file sample. reo 



reduce < sample. rei > sample. reo 



Just like running own file you can read packages that implement new commands or whole 
calculi with the command load_package. Table 2.1 briefly displays the packages available 
for Reduce. In principle, one can think of packages as ordinary Reduce code which is 
precompiled into a fast loading file. (Such files can be generated with the commands 
f aslout "filename" ; f aslend; ). Later, we will describe the Excalc package implementing 
the exterior calculus. You find a list of all packages and links to online references in the 
table of contents of the Reduce user's manual at [6] 

Also, the ZIB (Berlin, Germany) offers a set of references (in PDF format) for all packages 
at [9]. 
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ATP TATT 


integration for functions involving roots (James H. Davenport) 


A T> ATT TA /T 

AKIN u m 


algebraic numbers (Ebcrhard Schrfer) 


A CCTCT 


useful utilities for various applications (Hubert Caprasse) 


AVFPTni? 


veciui aigeoia ^uaviu naipei j 




package for computational commutative algebra (Hans-Gert Graebe) 


p A Al A T 


calculations in celestial mechanics (John Fitch) 


OriAlN LxrL VAK 


transformation of variables in differential equations (G. oluk) 


POA/TPAPT 
OVjlVlr ALy 1 


condensing of expressions with polynomial side relations (Anthony C. Hcarn) 




package for solving overdetcrmincd systems of PDEs or ODEs (Andreas Brand, Thomas Wolf) 


C Vll 


Dirac gamma matrices (V.Ilyin, A.Kryukov, A.Rodionov, A.Taranov) 




differential ecjuations and singularities (C. Dicrescenzo, F. Richard-Jung, E. Tournier) 


UYP AT P 
I.A LA I A 


calculus for differential geometry (Eberhard Schrfer) 


TPTFlF 
r ILTj 


code generation for finite difference schemes (Richard Liska) 


vjUjIN J_ it A IN 


r-nr\& ircnorotiATi in FHRTR A AT R /'Rot-Koto f^a-faoA 

cocic generation m r ^xii xvA-In , xi-A-irvjrt, k^j ^oaruara v^atesj 


PATTTDT 


display of functions and surfaces (Herbert Melenk) 


T> /A T? TU AT T? T) 


computation in multivariate polynomial ideals (Herbert Melenk, H.Michael Mllcr, Winfricd Neun) 




high energy physics (Anthony C. Hearn) 


TTATT AT Q 


Arithmetic for polynomial ideals (Herbert Melenk) 




Laplace and inverse Laplace transform (C. Kazasov et al.) 


T TT? 


functions for the classification of real n-dimcnsional Lie algebras (Carsten, Franziska Schbcl) 


T TIVTTTQ 


a package for finding limits (Stanley L. Kamcny) 


IjIIN TIN Hi Vol 


linear inequalities and linear programming (Herbert Melenk) 


TsITTATTTRTP 
in u ivir^rtiv,' 


solving numerical problems using rounded mode (Herbert Melenk) 


\JLJFjO\J1jv Hi 


ordinary differential equations (Malcolm MacCallum et al.) 


( \ \y\~ \ 1 CWJ LY 1 
v 7 1 \ 1 riU V 


calculus for scalar and vector quantities (J.W. Eastwood) 


111 i DUr 


additional support for non-commuting quantities (Mathias Warns) 


.r 1V1 


general algebraic pattern matcher (Kevin Mclsaac) 


Dr APTPP1M 
iv I ^ A v l rljV^iN 


manipulation of chemical reaction systems (Herbert Melenk) 


~d T T-TT TPT 
XtJ-JX 1 1, 1 1\1 


TcX and LaTeX output (Richard Liska, Ladislav Drska, Werner Antwciler) 


1 A 1 \ ) 1 o 


roots of polynomials (Stanley L. Kamcny) 


>jv All I. 


optimization of numerical programs (J. A. van Hulzen) 


or Drj 


symmetry analysis for partial differential equations (Fritz Schwarz) 


QlDTPr^T-PAT 
orrj^r IN 


package for special functions (Chris Cannam et al.) 


Q "D J? T? AT O 


package for special special functions (Victor Adamchik, Winfricd Neun) 


SYMMETRY 


symmetry-adapted bases and block diagonal forms of symmetric matrices (Karin Gatermann) 


SUM 


sum and product of series (Fuji Kako) 


TAYLOR 


multivariate Taylor series (Rainer Schpf) 


TPS 


univariate Taylor series with indefinite order (Alan Barnes, Julian Padget) 


WU 


Wu Algorithm for polynomial systems (Russell Bradford) 



Table 2.1: Reduce packages 
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2.2 



Examples 



2.2.1 



L Generating a Julia set 

7, file "red3" 



on comp$ 



procedure julia(c,s,f ile)$ begin 
on complex, rounded$ 
precision(6)$ 
1:={}$ 

for x:=-3/2*s:3/2*s do « write x$ for y : =-3/2*s : 3/2*s do « 
z:=x/s+i*y/s$ 

j:=0$ repeat z:=z**2+c until 

abs(repart(z))>2 or abs (impart (z) ) >2 or (j:=j+l)=50$ 
if j=50 then l:={x,y}.l$ 
» »$ 

off complex, rounded$ 
out file$ 

for each point in 1 do 

write first (point) , " ", second (point) $ 
shut file$ 
end$ 

showtime$ 

julia(-0 . 11+0 . 67*i , 100 , "12 . red" ) $ 
showtime ; 

7, time: 780730 ms 

Vein gnuplot : 

7.set data style dots 

7.set nokey; set noxtics; set noytics; set size square; set noborder 

•/.plot "12. red" 

7.set terminal postscript 

7set output "julial.ps" 

7replot 

2.2.2 Proving the Casimirs of the Poincare group 

7. file "red2" 

7 defining the epsilon and the Minkowsi metric: 
operator eps$ antisymmetric eps$ 
let eps(0,l,2,3)=>l$ 
array mink (3, 3)$ 

mink(0,0) :=1$ mink(l,l) :=mink(2,2) :=mink(3,3) :=-!$ 
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Figure 2.1: The Julia set for z — > + (— 0.11+0. 67i). (File julial.ps). 



7, declaring the generators of the Poincare group: 
operator j,p$ noncom j,p$ antisymmetric j$ 

% the algebra of the Poincare group and general rules for the Lie bracket: 
% note that the order of these rules is very important ! 
infix lie$ 
let { 

( p(~a) lie p(~b) ) => 0, 

( p(~a) lie j(~b,~c)) => mink(a,b) *p(c) - mink(a, c) *p(b) , 
(j(~a,~b) lie p(~c) ) => - (p(c) lie j(a,b)), 

(j(~a,~b) lie j(~c,~d)) => - mink(a,c)*j (b,d) - mink(b,d)*j (a,c) 

+ mink(a,d)*j (b,c) + mink(b,c)*j (a,d) , 
(- ~x lie ~y) => - (x lie y) , ( ~x lie - ~y ) => - (x lie y) , 
(~x lie ~y) =>x*y-y*x 

}$ 

°/« rules for applying commutator to order products of operators 
l:=oplist:={p(0) )P (l),p(2),p(3),j(0,l),j(0,2),j(0,3),j(2,3),j(3,l),j(l,2)}$ 
for each x in 1 do 

for each y in (l:=rest 1) do let y*x=>x*y- (x lie y)$ 

°/o the momentum square (p2) , the Pauli-Lubanski (pi) , and its square (pl2) : 
operator p2,pl,pl2$ 

p2:=for a:=0:3 sum for b:=0:3 sum mink(a,b)*p(a)*p(b)$ 
for a:=0:3 do pl(a) := (-1/2)* 

for b:=0:3 sum for c:=0:3 sum for d:=0:3 sum eps(a,b,c,d)*j (b,c)*p(d)$ 
pl2:=for a:=0:3 sum for b:=0:3 sum mink(a,b)*pl(a)*pl(b)$ 



°/ the commutators of p2 and pl2 with all generators : 
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for each x in oplist do write " [ p2 , ",x," ] = ",p2 lie x$ 
for each x in oplist do write "[ pl2 , ",x," ] = ",pl2 lie x; 



2.3 Advanced structures 



2.3.1 Lists 



Lists are ordered sets. They are constructed via 

{ elements }; 

The meaning of the commands to manipulate list should become clear from the following 
example: 

li:={a,b,c,d>; 

first li; '/„-> a 

second li; %-> b 

third li; '/,-> c 

part (li,4) ; 7,-> d 

rest li; '/,-> {b,c,d} 

reverse li; °/ -> {d,c,b,a} 

length li; '/-> 4 

append (li,{e,f ,g}) ; '/„-> {a,b,c,d,e,f ,g} 

. li; °/„-> {0,a,b,c,d> 

. 1 . 2 . 3 . li; '/.-> {0,l,2,3,a,b,c,d} 



2.3.2 Loops and conditions 

The syntax for the loop and conditional commands are 

for name := start: stop [do I sum I product I collect I join] statement; 

for name := start step stepsize until stop [do I sum I product I collect I join] statement; 

for each name in list [do I sum I product I collect I join] statement; 

while boolean do statement; 

repeat statement until boolean; 

if boolean then statement; 

if boolean then statement else statement; 



The do action simple executes the statement in each iteration. With the sum, product 
actions, the hole for statement returns the sum or product of all statements in the iteration. 
Analogously, the collect action returns a list of all statements and the join action returns 
the union of all statements (these have to be lists in this case!). 



Examples: 
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for x:=l:5 do write x; % writes the numbers from 1 to 5 

for each x in {2,4,8} sum x/2; '/,-> 7 

for each x in {a,b,c> join {x,2*x}; '/„-> {a,2*a,b,2*b,c,2*c} 

if x neq and not x>0 then write "impossible!"; % returns nothing 

As we can see, a boolean is an expression formed by the logical operators = neq > >= <= 
< and or not. 



2.3.3 Groups 

Instead of executing only one statement in each loop or in some conditional case, we can 
execute an arbitrary set of statements by embracing them to a group. A simple way to 
form a group is 

« IsLstatement; 2ndstatement; ... ; last-statement » 

This group is again a statement with the value of the last statement (if one did not put a 
terminator after the last statement). A more robust way to form groups is 

begin [scalar local jnames;~\ statements [return value;] end; 

Here, locaLnames are such that have no effect outside of this group. They are initialized 
with 0. This group returns the specified value. To produce an output within a group one 
has to use the write command. 



2.3.4 Procedures 

Procedures are operators with an explicitly defined set of parameters and an explicitly 
defined functionality. The syntax for defining procedures is: 

procedure name (parameter_names) ; statement; 

Here, statement defines the (return) value of the procedure and is usually a group statement. 

2.3.5 Operators 

Operators can be declared prefix or infix by: 

operator names; 
infix names; 

The declaration 

precedence name, next_lower_precedence-operator; 
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specifies the precedence of this operator: The declared operator has just higher precedence 
than the operator specified. 

Operators can have one of the following properties: 

■ totally symmetric or antisymmetric, 

• even or odd with respect to parity of the parameter, 

■ linear, or 

• noncom-muting under the multiplication * 
All of these properties can be declared by 

property operator --names; 
where property is one of symmetric , antisymmetric, even, odd, linear, noncom. 



2.3.6 Arrays, matrices 

For handling multicomponent objects Reduce offers the following declarations: 

array ( dimensions- 1) ; 

matr ix(co-dimension, contra-dimensiori) ; 

Both declarations initialize all components to be zero. The array can have arbitrary many 
indices (slots) and they are counted starting at 0. The are no special operators defined for 
handling arrays - all operations have to be done by hand, mostly with for loops running 
over the indices. 

The matrix has exactly two indices which start counting at l! A matrix can also be 
constructed by the mat command where the components have to be structured in tuples: 
e.g. mat ( (1 , 2 , 3) , (4 , 5 , 6) ) ; constructs a 2 x 3-matrix that Reduce displays as 

[1 2 3] 
[ ] 
[4 5 6] 

Reduce offers the following routines to handle matrices: 

tp trace det rank mateigen (matrix) ; 

They calculate the transpose, trace, determinate, rank, and the eigenvalue equation and 
eigenvectors of the given matrix, respectively. Also, matrices can be added, multiplied, and 
divided. (A division means multiplication with the inverse.) A matrix can be constructed 
via the 
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2.4 Export / import 
2.4.1 Storing results 

Storing results is very important for large problems. Reduce offers no extra command for 
storing expressions but it is quite easy to write a little procedure that stores e.g. an arbitrary 
list of names together with their values. For this we write appropriate assignments into an 
extra file with the nat switch turned off. Such a procedure could read 

procedure store(f ilename,namelist,valuelist)$ begin 
off nat$ 
out filename$ 

for each name in namelist do << 

write name ,":=", first valuelist$ 

valuelist : =rest valuelist$ 
»$ 

write "end"$ 
shut filename$ 
on nat$ 
end$ 

With this procedure defined we can execute the lines 

f:=sin(x)$ g:=cos(x)$ 
store ( " FandG" , "f " , "g" , f , g) $ 
clear f ,g$ 
in "FandG"$ 

f;g; 

where, the filename "FandG", the name list {"f","g"}, and their values {f ,g} is passed 
to the procedure store. This produces the file "FandG" in the current directory reading 

f :=sin(x)$ 
g:=cos(x)$ 
end$ 

such that the in command reassigns their values to f and g. 



2.4.2 TeX 

To export algebraic expression to TeX you need to load the package tri. This packages 
provides two switches tex and texbreak which causes Reduce to format any output in 
TeX syntax (texbreak also breaks the line in long equations). Example: 

1 : load_package tri ; 

*** global 'metricu!*' cannot become fluid 
*** global ' indxl ! * ' cannot become fluid 
7, TeX-REDUCE-Interface 0.50 
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7« set greek asserted 

7« set lowercase asserted 

7o \tolerance 10 

7« \hsize= 150mm 

2: depend(f,x); 

3: df(f.x); 
$$ 

\f rac{\partial f}{ 
\partial x} 

$$ 

For some TeX-symbols Reduce uses own macros which can be included in TeX with an 
\input{tridef s} command, (You can find the tridefs.tex file at [5].) To export all the 
output of a Reduce file in one file in TeX format one could add the lines 

load_package tri$ 

off msg$ 

on texbreak$ 

out "sample . out .tex"$ 

at the beginning of the Reduce file and the following lines at the end 

shut "sample. out. tex"$ 
quit$ 

Turning off the msg switch prevents errors. For further information on the tri package see 
http : //www . uni-koeln . de/REDUCE/3 . 6/doc/ tri . ps 



2.4.3 Fortran 

Exporting to Fortran is very similar to exporting to TeX. Reduce offers the fort switch 
to display all output in Fortran syntax. Unfortunately, Reduce can not generate whole 
Fortran procedures - as Maple does, e.g.. 

2.4.4 Maple 

Reduce offers, of course, no commands to export expressions and equations to Maple. 
However, in order to profit form the merits of both systems it is advantageous to know 
how to transport expression between them. Here we discuss how to export from Reduce 
to Maple. 

The strategy I use is to collect all expression in one list. As first and last entry I insert 
" [null" and "null] ", which are the braces to form ordered lists in Maple. Then I intro- 
duce rules to replace expressions which Maple does not understand, e.g. let { f =>f (x) , 
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df (f ,x)=>dif f (f (x) ,x) };. After reformulating the list one switches nat off and writes 
the list into a file. See the following example to pass a list of two expressions to Maple 

depend(f ,x)$ 

el:=int(f ,x)+f+x**2$ 
e2:=df (f ,x)+sin(x*f/3)$ 

maplist : ={ " [null " , el , e2 , "null] " >$ 

operator f ,diff$ 

let f=f (x) ,df (f ,x)=diff (f (x) ,x)$ 
maplist : =maplist$ 

off nat$ 

out "sample .map" $ 

write "maplist : =op(" .maplist , ") ; "$ 
shut "sample .map" $ 
on nat ; 

which results in the file sample .map: 

maplist : =op ({ [null , 

f(x) + int(f(x),x) + x**2, 

diff(f(x),x) + sin((f (x)*x)/3) , 

null]»;$ 

This file can be read by Maple with the command read ("sample .map") :. Although the 
$ sign at the end of the file will cause an error message, the maplist is read correctly (the 
LaTeX code for the following line was produced by Maple with the latex command): 



To display functions with gnuplot and to save them as postscripts you need to load the 
package gnuplot. This package provides the commands plot , gnuplot , plotreset, and 
plot show. With plot you will automatically open a gnuplot window that displays your 
graph. Example: 

load_package gnuplot$ 

plot(cos(x)*cos(y) ,x=(-pi. .pi) ,y=(-pi. .pi) .contour) ; 

will open a window that looks like the one in figure 2.2 

The syntax of this plot command is similar to that in gnuplot. The command 
gnuplot (cmd,paraml ,param2, . . .) ; 




2.4.5 Gnuplot 
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Figure 2.2: The gnuplot window and the postscript produced by gnuplot 

is supposed to execute all other gnuplot commands (which it does not always!) The pa- 
rameters of the gnuplot command com have to be separated by commas (not by blanks as 
in gnuplot). Example: 

gnuplot (set , logscale , y) $ 
plot(x**2) ; 

will produce a logarithmic plot of x 2 . Finally plotting into a file in postscript format is 
achieved by adding the options terminal=postscript and output=/i/ename to the plot 
command. Example: 

plot(cos(x)*cos(y) ,x=(-pi. .pi) ,y=(-pi. .pi) , contour ,terminal=postscript ,output="gnutest . eps") ; 

produces the postscript displayed in figure 2.2. For more information on the gnuplot 

package see 

http : //www . zib . de/Symbolik/reduce/moredocs/gnuplot . pdf 
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Chapter 3 
Maple I 



3.1 The simplification principle 

Just as Reduce, Maple is an input-output machine. But, as main difference to Reduce, 
Maple will not automatically reformulate any input obeying some rules. Instead, Maple 
only executes explicit commands in the input but cites expressions without commands ver- 
bally. Hence, the simplification performance of Maple is controlled with explicit 
commands only - not with automatic reformulation rules. This nature of Maple 
might stem from its origin: Maple is programmed in C and the command syntax is strongly 
influenced by this language. There exist no rules or switches, Maple becomes active only 
when a command or an operator is called. 

Maple offers a variety of commands to reformulate expressions. Most important though 
is one command, simplify, which tries to apply many different simplification methods to 
the expression and can be controlled with some options. For example, Maple will simply 
cite an input (x**2-l)/(x+l) ; as 

x 2 -l 
x + 1 

but the answer on simplify ((x**2-l)/ (x+1)) ; is 

x-l 

The general syntax of simplify is 

simplify (expr [ , simplification-methods] [,assume= assumptions] ); 

Without specifying the simplification methods, simplify applies all methods. The pa- 
rameter simplification-methods can be the name of a user-defined simplification procedure 
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'simplify/name' (expr) or any of the build-in methods 

Bessell, BesselJ, BesselK, BesselY, D, Ei, GAMMA, RootOf, LambertW, dilog, exp, 
In, sqrt, polylog, pochhammer, trig (for trig functions) , hypergeom (for hypergeo- 
metrics) , radical (occurrence of fractional powers) , power (occurrence of powers) , exp, 
In, and atsign ("@" - for operators) 

One main problem with Maple's simplification principle is, that simplify can not handle 
unknown (e.g. package-defined or user-defined) objects as expressions. For a tensor, e.g., 
the simplification command has to be applied on each component separately. Here, Reduce 
is clearly more elegant. 



3.2 The interface 



When you start a graphical version of Maple ('xmaple') you are offered a worksheet dis- 
playing the first prompt > . Behind the prompt, Maple expects a command consisting out 
of an statement and a terminator (; or :). Again, the terminator decides whether Maple's 
response is displayed (;) or not (:). Note that Maple is case sensitive! All examples 
displayed in this chapter are collected in the file examples .mws at [5]. The LaTeX code to 
display the examples here was produced by Maple (see Import / export). 

Simplifying 

> (x**2-l)/(x+l) ; 

x 2 -l 



X + 1 

x — 1 



> simplify((x**2-l)/(x+l)) ; 

> expand((x+y)**3) ; 

x 3 + 3 x 2 y + 3 x y 2 + y 3 

> A : = [ (x~2-l) / (x+1) , 0] ; simplify (A) ; 



A := 



,x 2 - 1 



pO] 



x 

[x-1, 0] 

> A:=[0,0] ; simplify (A) ; 

A := [0, 0] 
[0, 0] 

> A:=array([[(x**2-l)/(x+l) ,0] , [0,0]]) ; simplify(A); 

c 2 - l 

A:= ' 



x + 1 


x-1 





> A:=array([[0,0] , [0,0]]) ; simplify(A) ; 
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x-l 


Note that Maple makes an error in the last line because it does not handle arrays correctly! 
The operator % has the value of the last response of Maple (cf. ws in Reduce), i.e. in our 
case Similar, °/„°/ and V/X represent the second and third last responses of Maple. 

Usually, after typing the terminator one presses ENTER to execute the command. How- 
ever, Maple offers to group commands in one so-called execution groups. To insert a 
line feed in this group one presses SHIFT+ENTER. Pressing ENTER somewhere in an 
execution group then executes all commands in this group. 



3.3 User definitions 
3.3.1 Names and assignments 

Assignments 

> my_namel23 : =int (f (x) ,x) ; 

my_namel23 := j f{x) dx 

> '$& "' *@#! ' :=int(g(x) ,x) ; 

$& "' *@#! := Jg{x)dx 

> my_namel23:='my_namel23' ; 

> '$& " ' *@#! ' := ' '$& "' *<§#! 1 ' ; 

my_namel23 := my_namel23 
$& "' * @#! : = $& "' * @#! 

> i:=0: j :=2:k:=l: A.i.j.k; 

A021 

As in Reduce any string (which is not reserved or a keyword, see table 3.1) may be used as 
a name. If you want some unusual string to be a name, the string has to be embraced by 
back(!)-quotes, e.g. '123 string' can be used as a name. The . sign allows to append 
numbers to a name. 

The assignment operator is name:=expr. As in Reduce, the lhs name is identified with the 
rhs expression. There are three methods to unassign a name x: 

- First with x: = 'x' ;. The embracing by single quotes(!) means an unevaluation of the 
enclosed string. Hence, x : = ' x ' ; means an assignment of the unevaluated string x back to 
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Pi I infinity gamma Catalan true false 

Digits Order constants libname printlevel lasterror 
status 

and by do done elif else end fi for from if in intersect 
local minus mod not od option options or proc quit read 
save stop then to union while 

Table 3.1: Reserved names in Maple, 
x. 

- Second with unassign( 'x' ) ;. The unassign command also needs an unevalued string 
as parameter. 

- Third with restart; which, however, unassigns all names. This methods is most secure 
to free the memory resource maple allocated. 

As in Reduce it is important to distinguish between the expression assigned to a name and 
the value that Maple displays when evaluating the name. The following example shows 
the importance of the reassignment a:=a:. See also section 1.3.1. 

> a:=b: a; 

b 

> b:=l: a; 

1 

> a:=a:b:=2: a; 

1 

> a:='a'; 

> b:='b'; 

a := a 
b := b 

3.3.2 Aliases 

Aliases 

> alias(I='I') :macro(I='I') : 2*1; sqrt(-4); 

21 
2v^T 

> macro (I=sqrt (-1) ) : 2*1; sqrt(-4); 

2v^T 
2v^T 

> macro (1= 'I') : 2*1; sqrt(-4); 

21 
2v^T 



constants 
internal variables 

keywords 
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> alias (I=sqrt (-1) ) : 2*1; sqrt(-4); 

21 
21 

> alias(I='I') : 2*1; sqrt(-4); 

21 
2 V^T 

Maple has two possibilities to define aliases which are quite useful (in contrast to the 
define command in Reduce). The first is macro and the second is alias. They have the 
same syntax but they have different meaning: 

The macro (name=expr) ; command does not evaluate the expression expr before assigning 
it to the macro's name. Whenever this name will appear in a future expression the first 
thing Maple does is to replace it by the associated value expr. However, Maple will not 
replace expr by the macro name in outputs. The point of a macro is to introduce shorthand 
notations for the input. 

The alias (name= expr) ; command evaluates the expression expr before assigning it to the 
name. The value of expr cannot be a number. The effect of an alias is first, that Maple 
replaces all aliases by its value before evaluating expressions (as it does with macros), but 
second, that Maple also replaces expressions in the output by its alias if these exists one. 
So, the main point of an alias is introduce abbreviations in the input and output. The best 
example: The complex unit I is defined in Maple via: alias (I=sqrt (-1) ) ; 

Macros and aliases can be changed and can be unassigned by assigning them back to their 
names: macro (name=name) ; and alias (name=name) ; 



3.3.3 Substitutions and assumptions 

Substitutions 

> f:=cos(x): x:=0: f; x:='x' : 

1 

> f:=cos(x): simplify (subs (x=0,f )) ; 

1 

> f: = 'f: 

As in Reduce, one way to substitute identities in expressions is to make a local definition 
that represents the identity and let Maple reevaluate the expression. After this, f is still 
assigned to cos(x) but was evaluated to be 1. However, one should clear the assignment 
for x again with x : = ' x ' : . A more elegant way to substitute is the subs command the 
syntax of which is: 

subs (eqnjist, expr) ; 
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Eventually, one has to simplify again to take the substitution into account. 

Maple has one quite surprising and useful tool that Reduce lacks. Maple allows to make 
assumptions on unassigned names! With this, e.g., Maple simplifies \Zx* to x if one 
assumed x > before: 

Assumptions 

> sqrt(x~2); 

Vx* 

> assume (x>0) : about (x); 
Originally x, renamed x~ : 

is assumed to be: RealRange (Open (0) , infinity) 

> sqrt(x~2); 

x~ 

> x := 'x' ; 

x := x 

> sqrt (x~2) ; 

The syntax is assume (properties) ; The twiddle ~ behind the name x indicates that there 
are properties assumed for x. You can ask for properties of names with about and is and 
add more properties with additionally. Properties are deleted with the reassignment 
x: = 'x' ; 



3.3.4 Operators 

Operators 

> def ine(log2,log2(2)=l) : 

> log2(2)*log2(3); 

log2(3) 

> def ine(f ac,f ac(0)=l,fac(n: :posint)=n*f ac(n-l) ) : 

> fac(5); 

120 

> unassign(lo g 2,f ac) : 

Maple offers the define commands to declare new operators. The general syntax is 

define (name , properties) ; 
definemore (name, properties) ; 
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Here, properties are equations (rules) or any of the following key properties: linear, 
multilinear, orderless, flat (means associative) . Note that it is possible to specify 
the type (see sections 4.2.2 and 4.2.3) of a parameter (e.g. n: :posint to allow only positive 
integers for the parameter n). The command definemore adds more properties to an 
already defined operator. 

3.3.5 Functionals 

> f(x):=cos(x): f(x); f(0); f(y); 

cos(x) 
f(0) 
%) 

> f (x) : = 'f (x) ' : 
Functionals 

> g := x -> cos(x); g(x) ; g(0) ; g(y) ; D(g); 

g := cos 
cos(cc) 
1 

cos(y) 
—sin 

> h := (x,y) -> sin(x)*sin(y); D[l](h); D [2] (h) ; 

h := (x, y) — > sin(x) sin(y) 
(x, y) -> cos(x) sin(y) 
(x, y) -> sin(x) cos(y) 

> k := h (x -> (x,exp(x))); k(x) ; 

k := h@(x -» (x, e x )) 
sin(x) sin(e x ) 

> sinOexp; (sinOexp) (x) ; D(sinOexp); (D(sin@exp)) (x) ; 

sin@exp 
sin{e x ) 
cos@exp exp 
cos(e x ) e x 

> sin@@3; (sin@@3) (x) ; 

sin' 3 ) 
(sin( 3 ))(x) 

> g: = 'g': f: = 'f: h: = 'h' : k: = 'k' : 
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3.4 Solving 

Maple offers the solve command for solving equations. The syntax is 
solve (eqnJist, nameJist) ; 

For solving differential equations, use dsolve. For purely floating-point solutions use 
f solve. Use isolve to solve for integer solutions, msolve to solve modulo a prime; rsolve 
for recurrences, and linalg [linsolve] to solve matrix equations. 

Solving 

> solve (x~4-5*x~2+6*x=2, x) ; 

-1 + V3, -1 - Vs, 1, 1 

> dif f (y (x) ,x) -y (x) ~2+y (x) *sin(x) -cos (x) ; dsolve(7 n ); 

_ y( x ) 2 + y( x ) sin ( x ) _ cos ( x ) 

e (-cos(a:)) 

y(x) = -j h sin(a;) 

_Ci + / e (-cos(x)) ,j,x 

> diff (x(t) ,t)=y(t) ; diff (y(t),t)=x(t)+y(t); x(0)=2,y(0)=l; 

ix(i)=y(t) 
&y(t)=x(t)+y(t) 
x(0) = 2, y(0) = 1 

> res := dsolve ({7 7 7 , '/'/<,, 7o}, {x(t) ,y(t)}, type=mimeric , output=procedurelist) ; 

res := proc(rkf45 .x) . . . end 

> res(l); 

[t = 1, x(t) = 5.582168689244844, y(t) = 7.826891137110794] 
We see how to plot solutions (esp. of dsolve) in chapter 4. 

3.5 Commands and references 

Maple has many, many command - a lot more than Reduce. Here, I only tried to collect 
the most basic ones in table 3.2. Other than Reduce, Maple offers complete commands for 
almost any standard problem instead of providing only elementary commands that suffice 
to write procedures. Many of them might provide quick help. Unfortunately, this means a 
vast number of commands which can hardly be overviewed. Most of these one will never 
use anyway because they are too specialized. 

Because of these innumerable commands Maple offers the most efficient way to work with 
Maple is by using the build-in help browser. This browser is opened via the programs 
menu or, more efficiently, by positioning the cursor on a word and pressing 'Crtl+Fl'. To 
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begin with, it is useful to have a look on the 'Index...' topic (you find in the upper left 
window). There you find lists of expressions, functions, packages, etc. Alternatively, the 
'Topic search...' in the 'Help' menu it useful to search for commands etc. alphabetically. 

There is also online help in the internet available at [8] . This page is the producers address. 
You can find more addresses at [5]. 

There is also help available with in the non-graphical interface. The command ^subject 
will display help and examples for the specified subject. Also the commands index info 
usage example related {subject) ; provide help on the subject. 

Integration / differentiation /dependencies: 

> diff (cos(5*x) ,x) ; 

> int (sin(x) **2 ,x) ; 

—5 sin(5 x) 

1 / \ * / \ 1 
— - cos(x) sin(x) + -i 



> diff (f (x) ,x) ; diff(f,x); int (f (x) ,x) ; int(f,x); 

d_ t 

dx 



^cos(x)(x) 
sin(a;) 

cos(x)(:r) dx 

sin(x) 



/ 



> f(x):=cos(x): f(x); f(0); f(y); 

cos(x) 
f(0) 
%) 

> f (x) : = 'f (x) ' : 

Procedures 

> my_f actorial : =proc (k) 

> local x,l: x:=l: 

> for 1 from 1 to k do x:=x*l: od: 

> end; 

my .factorial := proc(fc) local 1 ; for (tofcdoi:=sx (od end 

> my_f actorial (16) ; 

20922789888000 

> unas sign (my _f actorial) : 

Selection, sets, lists 

> op( int(f(x),x) ); op( [1,2,3] ); 

f(x), x 
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1, 2, 3 

> lhs(x=l); rhs(f (x)=sin(x)) ; 

x 

sin(x) 

>{2,1}; [2,1]; 

{1,2} 
[2, 1] 

> seq( x[i], i=1..5 ); x[i] $ i=1..5; 

Xi, x 2 , x 3 , Xi, x 5 
Xi, x 2 , x 3 , Xi, x 5 
>{1,2} union {3,4}; {1,2} intersect {1,3}; 

{1, 2, 3, 4} 
{1} 

> f := x -> x**2: 

> map(f , [x,y] ) ; map(f ,x+(y*z) ) ; map(f ,x* (y+z) ) ; 



x 2 (y + z) 2 

Logical operators 

> (x>3) and (y<=l) ; 

3 - x < and y - 1 < 

> subs({x=5,y=0},7„) ; 

-2 < and - 1 < 

> eval ('/.); 

true 

Inert operators 

> Sum(l/2"i,i=l. .k) = sum(l/2"i,i=l. .k) ; 

E^- 2 (^ +1) + i 

> Limit (Sum(l/2~i , i=l . .k) ,k=infinity) = 

> limit (sum(l/2~i , i=l . .k) ,k=inf inity) ; 

k 

um y _L = i 

> value (Limit (Sum ( 1/2" i , i=l . .k) ,k=inf inity) ) ; 

1 

Float evaluation 

> evalf (Pi, 70) ; 

3.141592653589793238462643383279502884197169399375105820974944592307816 
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Matrices 

> with(linalg) : A : =matrix( [ [1 , 2 ,4] , [2, 1 ,2] , [4,2, 1] ] ) ; 
Warning, new definition for norm 

Warning, new definition for trace 

A :-- 

> det(A); 

> inverse (A); 



1 2 4 

2 1 2 
4 2 1 



r -1 


2 



2 


~3~ 

2 


3 

-5 


3 



2 


3 

-1 


3 


IT - 
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system control 

assignment 

unassignment 

substitution 

elementary 

logical 

composition 

operators 

selection 

sets/lists 

arrays 

loops 

conditions 
procedures 

input output 



plotting 

evaluation 

solving 

handing types 



restart % °/„°/ °/ °/„°/„ quit; # comment; with {package); 
readlib (command_name) ; interface print [f] 

name := expr; alias macro (name=expr) ; assume (properties); 

unassign (name); x: = 'x'; 

subs asubs (eqnjist, expr) ; subsop (eqn_list,expr) ; 

expr +-*/** expr; factorial (integer) ; expr mod integer; 

expr <<=<> = >>= expr; boolean and not or boolean; 

functional @ functional; functional @@ power; 

define (name , properties) ; 

op ([operand-number ,~\ expr) ; denom numer (expr); 
lhs rhs (eqn) ; collect [t]coeff[s] 

{ set^elements }; [ list_elements ]; seq (expr , range) ; 
expr $ range; set intersect minus union set; 
member map select 

array (ranges , component Jist) ; evalm (matrix-expression); &* 
&*() 

(for from by [to I while] do od) (for in do od) (while do od) 
break next 

(if then elif else fi) 

(proc local global options end) 

ERROR WARNING (message) ; RETURN (exprs) ; 

latex fortran C (expr); writeto appendto ("filename"); 
save names /'filename" ; read "filename"; 

readdata( "filename" [, format] ) ; writedata ("filename" Jist) ; 
fopen fprintf fscanf fclose 

plot display odeplot textplot plotsetup interface 
E/eval eval [a| gf |m | b | c | hf | f | r] value allvalues expr; 
solve [i |m| r | f | d] solve ( eqn Jist, name Jist) ; testeq (eqn); 
whattype (expr); type convert (expr , type) ; 



Table 3.2: Elementary commands in Maple. We do not display the syntax of all commands 
and the innumerous functional operators sin etc. For explanations see table 1.2. The 
slash in E/eval means that there exist both commands Eval and eval, the first of which 
is the inert version (non-executing) whereas the second tries to execute the operation. The 
construction eval [a I gf |m| . .] means that there exists evala evalgf evalm etc. 
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4.1 Examples 



4.1.1 Generating a Julia set 

# file " julia.map" 
restart : 

julia:=proc(c,s,f ile) 
local x,y,z,i,l: 
1:=[] :z:=0: 

for x from -3/2*s to 3/2*s do print (x): for y from -3/2*s to 3/2*s do 
z : =x/s+I*y/s : 

for i from while i<100 and abs(Re(z))<2 and abs(Im(z))<2 do z:=z**2+c: od: 
if i=100 then 1 : = [op(l) , [x,y] ] : fi: 
od: od: 

save l.file: 
end: 

time (julia(0. 3-0. 5*1, 200, "16. map")) ; 

with(plots) : 
read "16. map" : 

plot setup (ps ,plotoutput=" julia3 .ps" , plot opt ions="noborder") : 
listplot (1 , style=P0INT , axes=N0NE , symbol=PDINT , scaling=CDNSTRAINED) ; 
plotsetup (default) : 
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Figure 4.1: The Julia set for z — > + (0.3-0.5i). (File julia3.ps). 

4.2 Advanced structures 

4.2.1 Loops and conditions 

The syntax for loop and conditional commands are 

if boolean then statements [elif boolean then statements'] [else statements'] fi; 

for name from siart [by siep] to stop do statements od; 

for name from start [by step] while boolean do statements od; 

for name in Zis£ do statements od; 

while boolean do statements od; 

4.2.2 Procedures 

The syntax for defining procedures in Maple is 

name := proc (parameter si: -.types]) local names; global names; statements end; 
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The last of the statements needs no terminator and its value is automatically the return 
value of the procedure. One can specify the type of each parameter by adding : : type to 
the parameters name. Example: 

f := proc(n: : integer) factorial (n) end; 

only excepts integers as parameters. Error messages or warnings can be produced with 
ERROR and WARNING. To use other procedures within a procedure it is sometimes useful to 
catch their error messages to prevent the procedure to be aborted. This is done be using 
the traperror command. 

4.2.3 Types 

There are innumerous different types defined in Maple. The main point of these is that 

procedures and operators (also user defined ones) can check if the given parameters have 

the correct type and throw error messages if not. E.g., if you pass an integer to the display 

command display (5) ; it says: 

Error, (in display) invalid argument, 5, must be a plot structure, or a list/set/array 
thereof 

This error message if quite useful. Also, if you want to tinker plot structures yourself 
and pass them to display, it is of course useful to find the exact description/syntax of the 
plot type in Maple's help page plot, structure. With this knowledge you could also 
manipulate plots. 

The commands 

whattype ( expr) ; 
type ( expr, type) ; 
convert ( expr, type) ; 

query the type of expr, check if expr has the type type, and convert expr such as to have 
the type type, respectively. See the help for the type for a list of all types Maple provides. 

4.3 Export / import 

For export and import there are four groups of commands. The first is very useful to store 

expressions. They are: 

save names," filename" ; 
read "filename" ; 

Example: 

f :=dif f (sin(x) ,x) ; 
save f , "sample .mai" ; 
restart ; 

read "sample .mai" ; 
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Here, the save command produces the file sample. mai which simply reads: 
f := cos(x) ; 

The following restart clears all assignments and cleans the memory buffer. The read 
command read the file as if it were typed in. 

The second group of commands is similar to the out command in Reduce and only write 
into files. The commands are 

writeto ("filename") ; 
appendto( "/Rename") ; 
writeto (terminal) ; 

They send all Maple's output to the specified file, writeto (terminal) ; closes the file and 
redirects the output to the worksheet. 

The third group of commands allows to read and write ordinary ascii data sheets as. Lists 
are used as buffer of such files: 

readdata( "/i/ename" [, format-list] [, size-of-one-data-tuple]) ; 
writedata( "filename" , list) ; 

readdata returns the data as list. The format-list can, e.g., be [integer , float , integer] 
for producing 3-tuples of data. 

The fourth group of commands implements the syntax of the programming language C. 
They are 

f := f open ("filename", READ I WRITE I APPEND, TEXT I BINARY) ; 
f printf (f , structurestring , exprs) ; 
f scanf (f , structurestring, variables) ; 
f close(f ) ; 

The meaning should be clear if one knows C. The commands feof , f flush, f remove, 
filepos also belong to this group. With these commands you can produce files with 
arbitrary content. 

4.3.1 LaTeX, Fortran, and C 

The commands latex (expr) , f ortran(eicpr) , and C(expr) are nice facilities to export to 
LaTeX, Fortran of the programming language C. See these examples: 

> restart: readlib(C) : 

> latex(Limit (Sum(l/2~i , i=l . .k) ,k=inf inity) ) ; 

lim _{k\rightarrow \infty }\sum _{i=l}~{k}\lef t ({2}~{i}\right )~{-l} 

> f := x -> Pi*ln(x~2)-sqrt(2)*ln(x~2)~2; 
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> C(f .optimized) ; 

/* The options were : operatorarrow */ 

#include <math.h> 

double f (x) 

double x; 

{ 

double t2; 
double tl; 
double t4; 
double t5; 
{ 

tl = x*x; 
t2 = log(tl) ; 
t4 = sqrt(2.0) ; 
t5 = t2*t2; 

return(0.3141592653589793El*t2-t4*t5) ; 

} 

} 

> f astgrow:=proc(k) 

> local x,i; 

> x:=l; 

> for i from 1 to k do x:=x+x**x; od; 

> x; 

> end; 

> f ortran(f astgrow) ; 

real function fastgrow(k) 
real k 

integer i 
real x 

x = 1 

do 1000 i = l,k,l 
X = x+x**x 
1000 continue 

f astgrow = x 
return 
end 



A possibility to export whole Maple worksheet to Latex is with the menu File - Export 

As - LaTeX. I used this method in this paper. This will produce a LaTeX document 

reproducing the appearance of the worksheet. The Maple style files have to be installed 

for this. One can get them at 

http : / /www . maplesof t . com/latex_patch . html 

and provisionally copy them in the same directory as the LaTeX file. 
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Figure 4.2: Plot produced by Maple 



4.3.2 Plotting 



Again, Maple offers innumerous commands to generate plots. Here, I present only those 
that I find most elementary: 

plot (functionalist [ , horizontaLrange] [,verticaLrange] [ , plot_options] ) ; 
textplot ({ [x-coord, y-coord, " text"] , . . } [ , plot-options] ) ; 
display ( plot-structure Jist [ , plot-options'] ) ; 
with(plots) ; 

odeplot (numeric-dsolve-output, function J,ist,xjrange [ , plot-options] ) ; 
plotsetup(c?ewce_ij/pe [.options] ) ; 

For handling plots it is important to realize that all plot commands return a plot structure 
of type PLOT. See this example: 

plt:=plot(sin(x) ,x=-Pi. .Pi) : 
save pit, "sample. pit" : 
restart : 

read " sample. pit " : 
display(plt) ; 

Here, we assign the plot of sin(x) to the variable pit and store this variable in the 
sample. pit file. This is very useful when it takes long to produce plots (e.g. with 
odeplot!). The following commands produce the postscript in figure 4.2 and write it 
in the file mapplot .ps. 

with(plots) : 

solution := dsolve({ diff(y(x),x) = sin(x*y(x)) ,y(0)=2},y(x) ,type=numeric) : 

pit := odeplot (solution, [x,y(x)] ,0. .6,labels=[x,y] ,labelfont= [TIMES, ITALIC, 22] ) : 

txt := textplot ( [2.5,2.5, "Solution of y'=sin(x*y), y(0)=2"] , align=RIGHT, font= [TIMES, ROMAN, 22] ) : 

plot setup (ps , plot output= "mapplot .ps" , plot opt ions="noborder") : 

display ({pit , txt} , axesf ont= [TIMES , ROMAN , 18] ) ; 

plotsetup (default) : 
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Tensor type calculi on Reduce and 
Maple 

5.1 Maple 

5.1.1 The tensor package 

Maple provides a the tensor package to implement the tensor calculus. This package is 
strongly oriented to perform the standard calculations in Einsteinian gravity, i.e. to calcu- 
late the Christoffel symbol, the Riemannian curvature tensor, the Ricci tensor, etc. from 
a given symmetric metric. For this special purpose the package is quite useful. Of course, 
the package provides also elementary operations to work with tensors: adding tensors 
lin_com, the tensor product prod, raising, lowering, contracting, and (anti-) symmetrizing 
indices raise, lower, contract, antisymmetrize, symmetrize. However, as one can 
see from table 5.1, the syntax of these commands is not at all in analogy to the usual 
notation of tensor calculus: To setup a general 3rd rank tensor A in 2D space with two 
covariant and one contravariant indices one has to declare 

restart : with(tensor) : 

A: =create ([-1,-1,1] , array ( 

[[[a.l.l.l(x.l,x.2),a.l.l.2(x.l,x.2)] , [a. 1 . 2 . 1 (x . 1 ,x . 2) , a. 1 . 2 . 2(x. 1 ,x. 2)] ] , 
[[a.2.1.1(x.l,x.2),a.2.1.2(x.l,x.2)] , [a. 2 . 2 . 1 (x . 1 ,x . 2) , a. 2 . 2 . 2(x. 1 ,x. 2)] ] ] 

)); 

To specify the tensor type, ' -1 ' means a covariant index whereas ' 1 ' means a contravariant 
(in contrast to what is said in the help page of create). Maple starts counting indices from 
1. As you see, one has to introduce general names for all components oneself (that include 
the dependency on coordinates). The following examples show other tensor manipulations: 
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declarations 
elementary ops 



differentiation 

more ops 

certain tensors 
coord trans 
GR tensors 

Newmann-P. 
extract eqns 



create (tensor_type, components-array) ; entermetric 

lin_com (coefficiants-and-tensors) ; 

prod ( lst-tensor, 2ndJ,ensor, index ..pairs -to '-contract) ; 

raise lower (metric, tens or, indices) ; 

contract (tensor , index-pairs) ; 

symmetrize ant i symmetrize (tensor , index_array) ; 
get_compts get_char get_rank (tensor) ; 

partial_diff exterior_dif f (tensor , coord-array) ; 
cov_dif f (tensor, coord-array, ChristoffeLsymbol) ; 
Lie_diff directional_dif f (tensor , vector , coord-array) ; 

compare invert dual permute_indices exterior_prod act 
conj frame commutator 

Levi_Civita dlmetric d2metric 

Jacobian change_basis transform 

tensorsGR display_allGR displayGR Christof f ell 
Christoffel2 Einstein Ricci Ricciscalar Riemann Weyl 
connexF RiemannF invars petrov 

convertNP npcurve npspin 

Killing_eqns geodesic_eqns 



Table 5.1: The commands of Maple's tensor package. The syntax is only specified for 
elementary commands. 



usual notation 


tensor package syntax 


5 Aij k + Bij k + 2 Cij k 


lin_com(5,A,B,2,c) ; 


A R * 

^ijk ±J lm 


prod(A,B, [1,3]); 




contract (A, [1,3] ) ; 


di Ajk 


partial_diff (A, [x. l,x.2,x.3,x.4] ) ; 



In contrast, we will see the Excalc allows tensor manipulations that are a direct mirror 
of the usual notation. The entermetric command allows to enter a metric interactively. 
Personally, I find this unacceptable for a CA system. 

Finally, we show how simple it is to calculate the standard GR tensors with the tensor 
package for a given metric. In the example we derive the Schwarzschild metric as solution 
of the Einstein equation for a general static, spherically symmetric ansatz. 



> restart : with (tensor) : 
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> gdd:=create( [-1,-1] , array (symmetric , [ 

> [f (r), 0,0,0] , 

> [0,-g(r)/f (r),0,0] , 

> [0,0,-r~2,0] , 

> [0,0,0,-r~2*sin(theta)~2] 

> ])); 

gdd := tablc([ 



compts 



f(r) 









S(r) 
f(r) 







-r 2 



-r 2 sin(6») 2 
index .char = [—1, —1] 
]) 

> tensorsGR( [t ,r ,theta,phi] ,gdd, 

> guu,gdet,chrisl,chris2,riem,ric,ric_scalar,einstein,weyl) : 

> eqs : =get_compts (ric) ; 

1 f(r) (-2 r g(r) f(r)) + r (& g(r)) (& f(r)) - 4g(r) (& f(r))) 



egs := 



0, - 



l(r) 2 r 



,0,0,0 



1 -2rg(r) f(r)) + r (£ g(r)) (& f(r)) + 4 ( « g(r)) f(r) - 4g(r) (£ f(r)) 



0,0, 



rf(r)g(r) 

f(r)) r g(r) - r g(r)) f(r) - 2 g(r) 2 + 2 f(r) g(r) 
2 g(r) 2 



,0 



0,0,0, 



1 sin(g) 2 (2 (J; f(r)) r g(r) - r g(r)) f(r) - 2 g(r) 2 + 2 f(r) g(r)) 



2 g(r) 2 
> dsolve({eqs[l,l] ,eqs [2, 2] } , {f (r) ,g(r)}) ; 

{i(r)=_C2 + --^},{g(r)=-Cl} 



5.1.2 The GRTensor package 
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> restart : 

> read "/home/mt/usr/maple/lib5/grii .m" : 

> grtensorO : 

> makeg(ss) : 

Makeg 2.0: GRTensor metric/basis entry utility 

To quit makeg, type 'exit' at any prompt . {\large } 

Do you wish to enter a 1) metric [g(dn,dn)] , 

2) line element [ds] , 

3) non-holonomic basis [e(l) . . .e(n)] , or 

4) NP tetrad [l,n,m,mbar] ? 

> 2: 

Enter coordinates as a LIST (e.g. [r ,theta,phi,t] ) : 

> [t ,r ,theta,phi] : 

Enter the line element using d [coord] to indicate differentials, 
(for example, r~2* (d [theta] "2 + sin(theta) ~2*d [phi] ~2) 
[Type 'exit' to quit makeg] 
ds~2 = 

> f(r)*d[t]~2 - g(r)/f (r)*d[r] "2 - r~2* (d [theta] "2 + sin(theta)"2 * 

> d[phi]~2): 



If there are any complex valued coordinates, constants or functions 
for this spacetime, please enter them as a SET ( eg. { z, psi } ). 



Complex quantities [def ault={>] : 

>{}: 

The values you have entered are 
Coordinates — [t, r, 9, (f>] 



9a b 



f(r) 







Metric : 


f(r; 









-r 









-r 2 sin(#) 5 
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You may choose to 0) Use the metric WITHOUT saving it, 

1) Save the metric as it is, 

2) Correct an element of the metric, 

3) Re-enter the metric, 

4) Add/ change constraint equations, 

5) Add a text description, or 

6) Abandon this metric and return to Maple . 



> 0: 



Calculated ds for ss (.020000 sec.) 



Default spacetime = ss 
For the ss spacetime : 
Coordinates 
x(itp) 

x a = [t, r, 9, <P] 
Line element 

ds 2 =i(r) dt 2 _ ^ r \ dr2 _ r 2 d Q2 _ r 2 sin (£)2 d( 

f(r) 



makegO completed. 

> grcalc( R(dn,dn) ): 



Calculated detg for ss (.010000 sec.) 
Calculated g(up,up) for ss (.020000 sec.) 
Calculated g(dn,dn,pdn) for ss (.020000 sec.) 
Calculated Chr (dn,dn,dn) for ss (.020000 sec.) 
Calculated Chr (dn,dn,up) for ss (.040000 sec.) 
Calculated R(dn,dn) for ss (.040000 sec.) 



CPU Time = .150 

> grdisplay( R(dn,dn) ) : 

For the ss spacetime : 
Covariant Ricci 
R(dn, dn) 

_ 1 f(r)((J:f(r))(|:g(r)) r - 2 f(r)) g(r) r - 4g(r) (£ f(r))) 
4 



R a b — 



,0,0,0 



°'4 



1 -2 f(r)) g(r) r + (£ f(r)) (|g(r)) r + 4 g(r)) f(r) - 4g(r) f(r)) 



0,0, 



f(r)g(r)r 

1 -2(|: f(r)) r g(r) + f(r) r (f g(r)) + 2 g(r) 2 - 2 f(r) g(r) 
2 g(r) 2 



,0,0 



0,0,0 



1 sin(g) 2 (-2 (jL f( r) ) r g(r ) + f( r ) r g(r )) + 2 g(r) 2 - 2 f(r) g(r)) 

2 g(r) 2 
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> eqs : ={grcomponent (R (dn.dn), [1 , 1] ) ,grcomponent (R (dn,dn), 



> [2,2])}: 

> dsolve(eqs,{f (r) ,g(r)}) ; 



{f(r) =_C2 + — }, {g(r) = _CJ} 



restart : 

read "/f iles/home/part3/mt/usr/maple/lib5/grii .m" : 
grtensor () : 
makeg(ss) : 
2: 

[t,r,theta,phi] : 

f(r)*d[t]~2 - g(r)/f (r)*d[r]~2 - r~2* (d [theta] "2 + sin(theta)~2 * d[phi]~2): 

O: 

0: 

grcalc( R(dn,dn) ): 
grdisplay( R(dn,dn) ): 

eqs : ={grcomponent (R (dn,dn), [1 , 1] ) ,grcomponent (R (dn,dn), [2,2])}: 
dsolve(eqs,{f (r) ,g(r)}) : 
latex ({'/„}); 



Who doubts about the beauty of the exterior calculus?! Reduce offers a package Excalc 
that provides very handy methods to handle exterior forms as well as tensors. Fortunately, 
this package implements elementary declarations and operators that let all freedom to the 
user to define quantities - it does not offer complete GR routines as the tensor of GRtensor 
packages in Maple that impose their notations and conventions. Excalc really generally 
implements the tensor and exterior calculus. Tensors can also be handled by Reduce 
without the Excalc package - you can find some explanation for this in the appendix. 
However, since with Excalc it is a lot easier and more beautiful, I advice to use the Excalc 
package to manipulate tensors. 

The Excalc package is loaded with loacLpackage excalc$. It provides the commands and 
operators summarized in table 5.2. The table should explain most commands. There is 
online help available at [7]. We now give more detailed explanations on the basic commands. 

The most important declaration is pf orm which allows to introduce exterior forms with 
additional tensor indices. E.g., with the following we introduce a 2rd rank tensor t a 
4th rank tensor r (r 1 ^ 1 ), and a 2-form p with one index (p l G A 2 ): 

pform t(i,j)=0,r(i,j ,k,l)=0,p(i)=2$ 
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declare forms/tensors 


pf orm name [ (tens or_indices)] = degree; 


declare vectors 


tvector names; 


declare dependence 


fdomain name=name( dependency) ; 


tensor symmetries 


index_symmetries tensor: [symmetric [in index-pairs]] 
[antisymmetric [in index-pairs]] ; 


get form degree 


exdegree (form) ; 


exterior derivative 


d (form) ; 


partial differentiation 


@ (f orm , coord) ; 


inner product 


vector _| form; 


hodge dual 


# (form) ; 


Lie derivative 


vector | _ form; 


variational derivative 


vardf (form,non_tensor_form) ; 


declare coframe 


coframe coframe_def [with metric metric-def 1 with 
signature signature] ; 


declare frame 


frame name; 


display coframe 


displayf rame ; 


indices 


indexrange numbers -orjnames; 


set space dimension 


spacedim number; 


tangent vector 


@ (coord) ; 


epsilon tensor 


eps (indices) ; 


Christoffel symbol 


riemannconx name; 


others 


noether keep nosum renosum noxpnd forder remforder 
xpnd 



Table 5.2: The commands of Excalc. A form is in general an exterior form declared with 
pf orm which can have tensor indices also. For the meaning of indexjpairs, coframe^def, 
metric-def and signature see the text. 



After this declaration, these tensors/forms can be refered to by X,(lst_indexjname,2ndjindexjname) , 
etc., i.e. we can use arbitrary names for the indices. Specifying the index range lets Excalc 
run the index names over this range in all expression. 

Putting a minus sign in front of the index name lowers the index, without a minus sign the 
index is a contravariant index. Excalc knows the summation convention summing 
over identical lower and upper tensor indices. E.g., the following examples calculate the 
trace U l and perform the assignment Uj := rkij k : 

1 : load_package excalc$ 
*** ~ redefined 

2: pform t(i, j)=0,r(i, j ,k,l)=0,p(i)=2$ 
3: t(-i,i); 



i 
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t 



4: t(i,j) :=r(-k,i,j ,k) ; 

i j i j k 

t := r 

k 

5 : indexrange , 1 , y , z$ 

6: t(-i,i); 

1 y z 

t + t + t + t 
1 y z 

7: t(0,0) :=a$ t(0,0) ; 

a 

9: t (-0,0); 


t 



If we have defined a metric Excalc automatically calculates the (contra-)covariant com- 
ponents from assignments to (co-)contravariant components. For a diagonal metric, the 
response to input 9 would then be g 00 a. Specifying the metric structure is a little bit 
involved. Roughly, the syntax to introduce a frame and metric is 

coframe coframe-def [with metric metric_def I with signature signature] ; 

The following examples should make clear what the parameters mean: 

7. file "excl" 
load_package excalc$ 

7« symbolic 3D coframe with Euclidean metric: 
coframe o (0) , o (1) , o (2) $ 

7 polar coframe with Euclidean metric: 
coframe o(r)=d r,o(p)=r * d phi$ 

7, Schwarzschild metric in anholonomic/diagonalized form: 
coframe o(0)=f * d t, 

o(l)=l/f * d r, 

o(2)=r * d theta, 

o(3)=r*sin(theta) * d phi 
with signature (+1,-1,-1,-1)$ 



7o symbolic 4D coframe with arbitrary metric: 
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operator gll,o$ °/« introduce general notations 

for i:=0:3 do for j:=i+l:3 do gll (i , j ) : =gll ( j , i) $ % symmetrize metric 
let { metric_def => for i:=0:3 sum for j:=0:3 sum gll(i, j)*o(i)*o(j) }; 
coframe o(0) ,o(l) ,o(2) ,o(3) with metric g = metric_def$ 

°/ holonomic 4D coframe with arbitrary metric and coordinates x(a) : 
operator gll,o$ 

for i:=0:3 do for j:=i+l:3 do gll (i , j ) : =gll ( j , i) $ 

let { metric_def => for i:=0:3 sum for j:=0:3 sum gll(i, j)*o(i)*o(j) }; 
coframe o(0)=d xO,o(l)=d xl,o(2)=d x2,o(3)=d x3 with metric g = metric_def$ 
pform x(i)=0$ °/« introduce coordinates as indexed 0-form defined as: 

x(-0):=xO$ x(-l):=xl$ x(-2):=x2$ x(-3):=x3$ 

frame e$ °/ generate a frame e(-a) dual to o(a) 

e(-a) _| o(-b); °/ -> gll(a,b) 

The first example introduces a coframe index by indices running over 0, 1, 2 in the Euclidean 
(because no metric is specified) 3D space. The second example introduces alphabetic(l) 
indices running over r, p and defines the anhonomic coframe with respect to holonomic 
coordinate derivatives (dr,d(f)). Thus, also r and phi are introduces as coordinates (e.g. 
enabling to write @ r for the radial vector). The 3rd example introduces the Schwarzschild 
geometry with diagonal metric and anholonomic coframe. 

The 4th example introduces a general coframe with arbitrary metric. The metric (low-low) 
components have been predefined as symmetric operator. Since the coframe command does 
not accept for constructions, the variable metric_def was used as alias for the metric 
expression. It is better to use let rules to predefine expressions that will be used in the 
coframe definition. The 5th example introduces a holonomic coframe with arbitrary metric, 
coordinate basis x(i), and holonomic frame @ x(a) which coincides with e(-a). 

The command 
frame e$ 

introduces the name e for the respective frame (defined via the metric by J cP = b~\ or 
e i \°j = 9ij)- Thus, the last example of the coframe definitions will yield 

e(-a) J o(-b) \7„ throws gll(a,b) 



Specifying the symmetries of tensors (or tensor indices on forms) saves Reduce from calcu- 
lating redundant components and thus saves time. Excalc provides the index_symmetries 
command for this. Note that if no metric is specified before the index_symmetries com- 
mand assumes an Euclidean metric! The syntax given in the table is 

index_symmetries tensor: [symmetric [in index-pairsll [antisymmetric [in index-pairs'] ] ; 

Here, indexjpairs can be {i,j} etc. but also pairs of pairs of indices, i.e. {{i , j} , {k,l}} 
which means that the pair {i,j} (anti-) commutes with the pair {k,l}. The following 
example will clarify this: 
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index_symmetries t(i,j): symmetric! 

index_symmetries r (k, 1 ,m,n) : symmetric in {k, 1} , {m,n} antisymmetric in {{k,l},{m,n}}$ 

These commands mean = P l and r klmn = r lkmn = r klnm = — r mnkl _ 

The rest of the commands works rather canonical and can be understood by looking at some 
examples: The first example implements the Schwarzschild metric in the tensor calculus: 

7. file "exc2" 

% The Schwarzschild metric in tensor formalism 
load_package excalc$ 

7, introducing the metric (the symbolic cof rame is of no importance) : 

depend (f ,r)$ depend (g,r)$ 

coframe o(0) ,o(l) ,o(2) ,o(3) with metric 

g = f**2*o(0)**2 -g**2/f **2*o(l)**2 -r**2*o(2)**2 -r**2*sin(theta)**2*o(3)**2$ 

°/ introducing coordinates with running indices: 
pform x(i)=0$ 

x(-0):=t$ x(-l):=r$ x(-2) :=theta$ x(-3):=phi$ 

% calculating the Christoffel symbol: 
pform chris(i, j ,k)=0$ 

index_symmetries chris(i, j ,k) : symmetric in {i,j}$ 
chris (-i , -j , -k) := 

(1/2)* (@(g(-j,-k),x(-i)) + <§(g(-k,-i),x(-j)) - <§(g(-i,-j),x(-k)))$ 

°/ the Riemannian curvature : 
pform riem(i, j ,k,l)=0$ 

index_symmetries riem(i, j ,k,l) : antisymmetric in {i,j},{k,l} 

symmetric in {{i, j},{k,l}}$ 

riem(-i,-j ,-k,l) := 

@(chris(-j ,-k,l) ,x(-i)) + chris (-i , -m, 1) *chris (-j , -k,m) 
- @(chris(-i,-k,l) ,x(-j)) - chris (-j , -m, 1) *chris (-i , -k,m) $ 

°/ the Ricci tensor: 
pform ric(i,j)=0$ 

index_symmetries ric(i,j): symmetric$ 
ric(-i,-j) := riem(-m, -i , -j ,m) $ 

% the Ricci scalar: 
pform ricscalar=0$ 
ricscalar := ric(-i.i); 

let {g=>l,f=>sqrt(l-mm/r)}$ 
ric(-i,-j) ; 

end$ 

The next example proves the Bianchi identity in a purely symbolic manner: 
°/„ file "exc3" 
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°/ Proving the Bianchi identity in the exterior calculus 

7. 

load_package excalc$ 
indexrange 0,1,2,3$ 

7« general introduction of a connection: 
pform gamma(a,b)=l$ 

7» definition of the curvature: 
pform curv(a,b)=2$ 

curv(-a,b) : =d gamma(-a,b) - gamma(-a, c) " gamma (-c ,b) $ 
7 Bianchi identity: 

d curv(-a,b) - gamma(-a,c) ~ curv(-c,b) + gamma(-c,b) " curv(-a.c); 
%-> yields zero ! 

end$ 

The next example implements the Taub-NUT solution with electric charge in the exterior 
calculus (in the Kaluza-Klein formalism): 

7, nut.exi, Marc Toussaint, Cologne 
7. 99-Dct-19 

7. 

7o The Taub-NUT solution with electric charge 

load_package excalc$ 
off gcd,exp$ 

/o COFRAME /////////////////////////////^//^/^///////////////////////^////////^^ 

clear t , r , thet a , phi , o , e , f , rho , dt au , ds igma , m , n , q$ 

pform f=0,rho=0$ 

fdomain f=f (r) ,rho=rho(r)$ 

let { dtau => d t - 2*n*cos (theta) *d phi, 
dsigma => (r**2+n**2) *d phi }$ 

dim: =5$ 

coframe o(0) = f/rho* dtau, 

o(l) = rho/f* d r, 

o(2) = rho * d theta, 

o(3) = 1/rho* sin (theta)* dsigma, 

o(5) = d w + r/rho**2* q*dtau 
with signature (1,-1,-1,-1,-1)$ 
frame e$ 



at aitciat7 °/ w w w °/ °/ w °/°/ °/°/ °/ ww °/v°/ w °/ v v w ww v v vvv v °/ vv / v °/v vv °/v°/ vv °/ v vvv vvvv °/ v vw v v 



f := sqrt(r**2-2*r*m-n**2+q**2/4)$ 
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rho := sqrt(r**2+n**2)$ 

7, field equation fflmfflfflfflfflnfflnfflfflfflfflifflfflnfflfflMmnffl 

7o Christoffel symbol: 

pform chrisl (a,b)=l$ 
Riemannconx chrisl$ 
chrisl(a,b) :=chrisl(b,a)$ 

7 Riemannian curvature and Einstein 3-form: 

pform rie2(a,b)=2 , einstein3(a)=dim-l$ 

rie2 (-a,b) : =d chrisl (-a, b) + chrisl (-c ,b) "chrisl (-a, c)$ 
einstein3(a) :=l/2*#(o(a) ~o(b) ~o(c)) ~ rie2(-b,-c)$ 

7, Result: 

einstein3 (a) ; '/,-> zero (except the (5 , 5) -compopnent) 

end; 



5.3 The tensor calculus in Reduce 

The tensor calculus is realized in Reduce by manipulating arrays with elementary com- 
mands. Usually, the first step is to introduce countable (i.e. indexed) coordinates as an 
operator (giving the coordinate for each index). Next one can introduce a metric and 
calculate the inverse to be able to raise and lower indices. Then tensors are introduced as 
arrays and all operations are performed by explicitly running over the indices (with the for 
command). Whether the array contains covariant or contravariant components is specified 
by adding the characters '1' for lower and V for upper indices to the name. The following 
example will follow this scheme to calculate the curvature of a general static, symmetric 
geometry in polar coordinates and thereby derive the Schwarzschild metric. 

7o indexed polar coordinates: 
operator x$ 

x(0):=t$ x(l):=r$ x(2):=theta$ x(3):=phi$ 

7« the metric (gll) with two lower and its inverse (guu) with two upper indices: 
array gll(3,3) ,guu(3,3)$ '/, arrays are initialized with zero 

7. static, spherically symmetric ansatz, signature: + - - - : 
depend (f,r)$ 
depend (g,r)$ 

gll(0,0) :=f **2$ gll(l,l):=- g**2/f**2$ gll (2 , 2) : =-r**2$ gll(3,3) : =-r**2* (sin(theta) ) **2$ 

7o calculating the inverse by using Reduce' s matrix calculus: 
matrix m(4,4)$ 

for i:=0:3 do for j:=i:3 do m(j+l,i+l) :=m(i+l, j+1) :=gll(i, j)$ 
m:=l/m$ 7, calculates the inverse 
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for i:=0:3 do for j:=i:3 do guu(j ,i) :=guu(i, j) :=m(i+l, j+l)$ 
clear m$ 

7, calculating the Christoffel symbol of type low-low-low and low-low-up: 
array chrislll(3,3,3) ,chrisllu(3,3,3)$ 
for i:=0:3 do for j:=i:3 do « 

for k:=0:3 do chrislll(j , i ,k) := chrislll(i , j ,k) := 

(1/2) * (df (gll(j,k),x(i)) + df (gll(k,i),x(j)) - df (gll(i,j),x(k)))$ 
for k:=0:3 do chrisllu( j , i ,k) := chrisllu(i , j ,k) := 
for m:=0:3 sum guu(k,m) * chrislll(i, j ,m)$ 

»$ 

'/, calculating the Riemannian curvature (neglecting symmetries -> inefficient!) 
array riemlllu(3,3,3,3)$ 

for i:=0:3 do for j:=0:3 do for k:=0:3 do for 1:=0:3 do 
riemlllu(i, j ,k,l) := 

df (chrisllu(j ,k,l) ,x(i)) - df (chrisllu(i ,k,l) ,x(j)) 
+ for m:=0:3 sum ( chrisllu(i ,m, 1) *chrisllu( j ,k,m) 

- chrisllu(j ,m,l)*chrisllu(i,k,m) )$ 

7« the Ricci tensor: 

array ricll(3,3)$ 

for i:=0:3 do for j:=0:3 do 

ricll(i,j) := for k:=0:3 sum riemlllu(k,i, j ,k)$ 

7« the Ricci scalar: 

ricscalar := for k:=0:3 sum for 1:=0:3 sum guu(l,k)*ricll(k,l)$ 
{ricll(0,0) ,ricll(l,l) ,ricll(2,2) ,ricll(3,3)> where {g=>l,f=>sqrt(l-mm/r)>; 
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